rm(list=ls())
library(ggplot2)
require(data.table)
library(dplyr)
library(tidyr)
require(parallel)
library(bdc)
library(taxadb)
library(traitdataform)
library(pbapply)
library(tidyverse)
library(readxl)
library(lme4)
library(coefplot)
library(sjPlot)
library(sjmisc)
library(effects)
library(rgdal)
library(maptools)
library(rgeos)
library(terra)
library(MuMIn)
library(rnaturalearthdata)
library(lsmeans)
library(GGally)
library(tidyterra)
# find accepted species names
source(here::here("R/clean_taxa_functions.r"))
# Malin's transformation for moving values slightly inward.
transform01 <- function(x) (x * (length(x) - 1) + 0.5) / (length(x))
# De Kort transformation
deKort_trans <- function(p){
p <- scale(p)
p <- (p - min(p, na.rm = T))/(max(p, na.rm = T)-min(p, na.rm = T))
return(p)
}
# Malin's transformation before a logit transformation
logit_trans <- function(p){
p <- transform01(p)
p <- log(p/(1-p))
return(p)
}
#function to select random effect
testSingAl=function(x,f,data,n=1:length(x)){
#function testing for random effect structure
#x=set of variables to test as random effect
#f=formula of the fixed structure of the model
#data= data base with all variables
#n= levels of variable combination to test (start from 1 when signle variable are tested to the size of the set of variables to test in combination)
b=1
if(is.null(x)==F){
for(i in n){
v1=combn(x,i)
if(i>1){
v2=apply(t(v1),1,paste,sep="",collapse="+")
}else{
v2=t(v1)[,1]
}
for(j in 1:length(v2)){
print(paste(i," : ",j,sep=""))
f1=paste(f,"+",v2[j],sep="")
lme1=lmer(as.formula(f1),data,REML=F,na.action="na.fail",control = lmerControl(optimizer ='optimx', optCtrl=list(method='nlminb')))
a=summary(lme1)$optinfo$conv$lme4$messages
if(is.null(a)==F){
a=paste(a,sep="",collapse="_")
}else{
a=NA
}
res=data.frame(test=v2[j],R2m=r.squaredGLMM(lme1)[[1]],R2c=r.squaredGLMM(lme1)[[2]],aic=AIC(lme1)[[1]],aicc=AICc(lme1)[[1]],singular=isSingular(lme1)[[1]],nV=i,source=grep("Source",v2[j])[1],warning=a[[1]])
if(b==1){
resOK=res
}else{
resOK=rbind(resOK,res)
}
b=b+1
}
}
}else{
f1=f
lme1=lmer(as.formula(f1),data,REML=F,na.action="na.fail",control = lmerControl(optimizer ='optimx', optCtrl=list(method='nlminb')))
a=summary(lme1)$optinfo$conv$lme4$messages
if(is.null(a)==F){
a=paste(a,sep="",collapse="_")
}else{
a=NA
}
resOK=data.frame(test=f,R2m=r.squaredGLMM(lme1)[[1]],R2c=r.squaredGLMM(lme1)[[2]],aic=AIC(lme1)[[1]],aicc=AICc(lme1)[[1]],singular=isSingular(lme1)[[1]],warning=a[[1]])
}
return(resOK)
#output:
#test= random structure
#R2m= marginal R2 of the model
#R2c= conditional R2 of the model
#aic= Akaike information criterion of the model
#aicc= Akaike information criterion with small-sample correction
#singular= output of the model singularity test (FALSE= no signularity; TRUE= singularity issue)
#nV= number of variables tested as random effect
#warning= inform for warning during the model fit (such as convergence issue)
}
splist <- read.csv(here::here("Data/splist.csv"), header = T)
# remove duplicated sp_names
splist <- splist %>%
dplyr::filter(!duplicated(scientificName))%>%
mutate(Kingdom=kingdom,
Phylum=phylum,
Class=class,
Order=order,
Family=family)%>%
dplyr::select(scientificName,Kingdom,Phylum,Class,Order,Family,db)
splist$Genus <- sapply(splist$scientificName, function(x){
strsplit(x," ")[[1]][1]
})
biov1 <- read.csv(here::here("Data/biov1_fixednames.csv"), header = T)
# Fix references in biov1
biov1$Article_ID <- gsub("[^0-9.-]", "", biov1$Article_ID)
biov1$Article_ID <- as.numeric(biov1$Article_ID)
biov1$sp_name_std_v1 <- gsub("_"," ",biov1$sp_name_std_v1)
biov1 <- biov1 %>%
dplyr::select(ID,Article_ID,Study_ID,
Type,Param,Trend,SHIFT,UNIT,DUR,
v.lat.mean,v.ele.mean,
START,END,Sampling,Uncertainty_Distribution,Uncertainty_Parameter,
N,Grain_size,Data,ID.area,
Phylum,Class,Order,Family,Genus,sp_name_std_v1,
group,ECO,Hemisphere) %>% # select columns
mutate(
Type = case_when(
Type=="HOR" ~ "LAT",
TRUE ~ as.character(Type)),
Data = case_when(
Data=="occurence-based" ~ "occurrence-based",
TRUE ~ as.character(Data)),
spp = sp_name_std_v1,
SHIFT_abs = abs(SHIFT),
velocity = ifelse(Type == "LAT", v.lat.mean, v.ele.mean),
vel_sign = ifelse(velocity>0,"pos","neg"))
biov1 <- biov1 %>%
dplyr::filter(!is.na(sp_name_std_v1))
all(biov1$sp_name_std_v1 %in% splist$scientificName)
## [1] TRUE
# from continuous to categorical
q1=quantile(biov1$START,probs=c(0,0.25,0.5,0.75,1))
biov1$StartF=cut(biov1$START,breaks=q1,include.lowest=T)
q1=quantile(biov1$ID.area,probs=c(0,0.25,0.5,0.75,1))
biov1$AreaF=cut(biov1$ID.area,breaks=q1,include.lowest=T)
q1=quantile(biov1$N,probs=c(0,0.25,0.5,0.75,1))
biov1$NtaxaF=cut(biov1$N,breaks=q1,include.lowest=T)
# add ID to obs
biov1$obs_ID <- paste0("S",1:nrow(biov1))
summary(biov1$velocity)
## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## -5.8328 0.5415 1.9642 2.2635 2.7543 14.5492 1450
summary(biov1$velocity[which(biov1$vel_sign=="pos")])
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.009396 0.541507 2.055833 2.411912 2.831932 14.549230
# Class species as Terrestrial, Marine or Freshwater
Terv1 <- unique(biov1$sp_name_std_v1[which(biov1$ECO == "T")])
Terrestrials = unique(c(Terv1))
Mar1 <- unique(biov1$sp_name_std_v1[which(biov1$ECO == "M")])
Marine = unique(c(Mar1))
# Freshwater fish
FFishv1 <- unique(biov1$sp_name_std_v1[(biov1$Class == "Actinopterygii" | biov1$Class == "Cephalaspidomorphi") & biov1$ECO == "T"])
FFish = unique(c(FFishv1))
# Marine fish
MFishv1 <- unique(biov1$sp_name_std_v1[biov1$Class == "Actinopterygii" & biov1$ECO == "M"])
MFish = unique(c(MFishv1))
splist$ECO = NA
splist$ECO[which(splist$scientificName %in% Terrestrials)] <- "T"
splist$ECO[which(splist$scientificName %in% Marine)] <- "M"
splist$ECO[which(splist$scientificName %in% MFish)] <- "M"
biov1$ECO = NA
biov1$ECO[which(biov1$sp_name_std_v1 %in% Terrestrials)] <- "T"
biov1$ECO[which(biov1$sp_name_std_v1 %in% Marine)] <- "M"
biov1$ECO[which(biov1$sp_name_std_v1 %in% MFish)] <- "M"
splist$Group = NA
splist$Group[which(splist$Class == "Phaeophyceae")] <- "Chromista"
splist$Kingdom[which(splist$Class == "Phaeophyceae")] <- "Chromista"
splist$Group[which(splist$Phylum == "Rhodophyta")] <- "Seaweed"
splist$Kingdom[which(splist$Phylum == "Rhodophyta")] <- "Plantae"
splist$Group[which(splist$Family == "Elminiidae")] <- "Barnacles"
splist$Group[which(splist$Kingdom == "Bacteria")] <- "Bacteria"
splist$Group[which(splist$Class == "Holothuroidea")] <- "Sea cucumber"
splist$Group[which(splist$Class == "Aves")] <- "Bird"
splist$Group[which(splist$Class == "Insecta")] <- "Insect"
splist$Group[which(splist$Class == "Mammalia")] <- "Mammal"
splist$Group[which(splist$Class == "Arachnida")] <- "Spider"
splist$Kingdom[which(splist$Kingdom == "Viridiplantae")] <- "Plantae"
splist$Kingdom[which(splist$Phylum == "Tracheophyta")] <- "Plantae"
splist$Group[which(splist$Kingdom == "Plantae")] <- "Plant"
splist$Group[which(splist$Class == "Hydrozoa")] <- "Hydrozoa"
splist$Group[which(splist$Class == "Anthozoa")] <- "Sea anemones and corals"
splist$Group[which(splist$Class == "Polychaeta")] <- "Polychaetes"
splist$Group[which(splist$Phylum == "Mollusca")] <- "Molluscs"
splist$Group[which(splist$Class == "Malacostraca")] <- "Crustacean"
splist$Group[which(splist$Class == "Hexanauplia")] <- "Crustacean"
splist$Group[which(splist$Class == "Maxillopoda")] <- "Crustacean"
splist$Group[which(splist$Class == "Ostracoda")] <- "Crustacean"
splist$Group[which(splist$Class == "Branchiopoda")] <- "Crustacean"
splist$Group[which(splist$Class == "Asteroidea")] <- "Starfish"
splist$Group[which(splist$Class == "Ascidiacea")] <- "Ascidians tunicates and sea squirts"
splist$Class[which(splist$Class == "Actinopteri")] <- "Actinopterygii"
splist$Group[which(splist$Class == "Actinopterygii")] <- "Fish"
splist$Group[which(splist$Class == "Elasmobranchii")] <- "Fish"
splist$Group[which(splist$Order == "Perciformes")] <- "Fish"
splist$Group[which(splist$Class == "Chondrichthyes")] <- "Fish"
splist$Group[which(splist$Class == "Holocephali")] <- "Fish"
splist$Group[which(splist$Class == "Cephalaspidomorphi")] <- "Fish"
splist$Group[which(splist$Class == "Echinoidea")] <- "Sea urchin"
splist$Group[which(splist$Class == "Crinoidea")] <- "Crinoid"
splist$Group[which(splist$Class == "Holothuroidea")] <- "Sea cucumber"
splist$Group[which(splist$Class == "Reptilia")] <- "Reptile"
splist$Group[which(splist$Order == "Squamata")] <- "Reptile"
splist$Group[which(splist$Class == "Ophiuroidea")] <- "Brittle stars"
splist$Group[which(splist$Class == "Chilopoda")] <- "Centipedes"
splist$Group[which(splist$Class == "Diplopoda")] <- "Millipedes"
splist$Group[which(splist$Class == "Amphibia")] <- "Amphibian"
splist$Group[which(splist$Kingdom == "Fungi")] <- "Fungi"
splist$Group[which(splist$Order == "Balanomorpha")] <- "Barnacles"
splist$Group[which(splist$Phylum == "Nematoda")] <- "Nematodes"
splist$Group[which(splist$Class == "Myxini")] <- "Hagfish"
splist$Group[which(splist$Kingdom == "Chromista")] <- "Chromista"
splist$Family[which(splist$Genus == "Dendrocopus")] <- "Picidae"
######################################
biov1 <- merge(biov1[,-which(names(biov1) %in% c("Phylum","Class","Order","Family","Genus","Group","ECO"))],
splist[,c("Kingdom","Phylum","Class","Order","Family","Genus","Group","ECO","scientificName")],
by.x = "sp_name_std_v1", by.y = "scientificName",
all.x = T)
Use only latitudinal and elevational shifts
# v1
biov1 <- biov1 %>%
dplyr::filter(Type %in% c("ELE","LAT")) # Use LAT ELE shifts
# splist
sps <- unique(biov1$sp_name)
splist <- splist %>% dplyr::filter(scientificName %in% sps)
if(any(grep("sp[.]",biov1$sp_reported_name_v1))){
biov1 <- biov1 %>% dplyr::filter(!grepl("sp[.]",sp_reported_name_v1))
}
if(any(grep("sp[.]",biov1$sp_name_std_v1))){
biov1 <- biov1 %>% dplyr::filter(!grepl("sp[.]",sp_name_std_v1))
}
if(any(grep("cf[.]",biov1$sp_reported_name_v1))){
biov1 <- biov1 %>% dplyr::filter(!grepl("cf[.]",sp_reported_name_v1))
}
if(any(grep("cf[.]",biov1$sp_name_std_v1))){
biov1 <- biov1 %>% dplyr::filter(!grepl("cf[.]",sp_name_std_v1))
}
#remove freshwater fishes
biov1 <- biov1[-which((biov1$Class == "Actinopterygii" | biov1$Group == "Fish") & biov1$ECO=="T"),]
#remove marine birds
biov1 <- biov1[-which(biov1$Class == "Aves" & biov1$ECO=="M"),]
if(any(grep("sp[.]",splist$scientificName))){
splist <- splist %>% dplyr::filter(!grepl("sp[.]",scientificName))
}
if(any(grep("sp[.]",splist$scientificName))){
splist <- splist %>% dplyr::filter(!grepl("sp[.]",scientificName))
}
if(any(grep("cf[.]",splist$scientificName))){
splist <- splist %>% dplyr::filter(!grepl("cf[.]",scientificName))
}
if(any(grep("cf[.]",splist$scientificName))){
splist <- splist %>% dplyr::filter(!grepl("cf[.]",scientificName))
}
unique(biov1$Data)
## [1] "abundance-based" "occurrence-based"
biov1$Data[biov1$Data=="occurence-based"] <- "occurrence-based"
biov1$Data <- factor(biov1$Data, levels = unique(biov1$Data))
table(biov1$Data)
##
## abundance-based occurrence-based
## 9496 20744
unique(biov1$Sampling)
## [1] "TWO" "CONTINUOUS" "MULTIPLE(continuous)"
## [4] "IRR"
biov1$Sampling = ifelse(biov1$Sampling %in% c("IRR","MULTIPLE(continuous)"),"MULTIPLE", biov1$Sampling)
biov1$Sampling <- ordered(biov1$Sampling,
levels = c("TWO","MULTIPLE","CONTINUOUS"))
table(biov1$Sampling)
##
## TWO MULTIPLE CONTINUOUS
## 25819 1251 3170
unique(biov1$Grain_size)
## [1] "small" "large" "moderate" "very_large" NA
biov1$Grain_size <- ifelse(biov1$Grain_size %in% c("large","very_large"),"large",biov1$Grain_size)
biov1$Grain_size <- ordered(biov1$Grain_size,
levels = c("small","moderate","large"))
table(biov1$Grain_size)
##
## small moderate large
## 9841 10971 9427
unique(biov1$Uncertainty_Distribution)
## [1] "RESAMPLING(same)" "RAW"
## [3] "RESAMPLING" "MODEL"
## [5] "RESAMPLING+MODEL" "MODEL+RESAMPLING(same)"
## [7] "RESAMPLING(same)+DETECTABILITY" "RESAMPLING(same)+MODEL"
## [9] "DETECTABILITY"
biov1$Uncertainty_Distribution <- ifelse(biov1$Uncertainty_Distribution %in% c("RESAMPLING","RESAMPLING(same)"),"RESAMPLING",
ifelse(biov1$Uncertainty_Distribution %in% c("MODEL","MODEL+RESAMPLING(same)","RESAMPLING+MODEL"),"MODEL",
ifelse(biov1$Uncertainty_Distribution %in% c("DETECTABILITY","RESAMPLING(same)+DETECTABILITY"),"DETECTABILITY",
biov1$Uncertainty_Distribution)))
biov1$Uncertainty_Distribution <- ifelse(biov1$Uncertainty_Distribution == "RAW","OPPORTUNISTIC","PROCESSED")
table(biov1$Uncertainty_Distribution)
##
## OPPORTUNISTIC PROCESSED
## 9532 20708
#transform study area
biov1$ID.area <- log(biov1$ID.area)
This will be used to merge genetic data with bioshifts based on distance of observations
# # The input file geodatabase
# fgdb <- "C:/Users/brunn/NextCloud/bioshifts_v1_v2/v1/Study_Areas_v1/Study_Areas.gdb"
#
# # List all feature classes in a file geodatabase
# fc_list <- ogrListLayers(fgdb)
#
# # Get centroids
# cl <- makeCluster(detectCores()-1)
# clusterExport(cl, c("readOGR", "gCentroid", "fgdb"))
#
# centroids <- pblapply(fc_list, function(x){
# fc <- readOGR(dsn=fgdb,layer=x)
# c <- data.frame(gCentroid(fc))
# names(c) <- c("centroid.x", "centroid.y")
#
# geomet <- data.frame(terra::geom(terra::vect(fc)))
#
# max_lat <- order(geomet[,"y"], decreasing = T)[1]
# max_lat <- geomet[max_lat,c("x","y")]
# names(max_lat) <- c("max_lat.x", "max_lat.y")
#
# min_lat <- order(geomet[,"y"])[1]
# min_lat <- data.frame(geomet[min_lat,c("x","y")])
# names(min_lat) <- c("min_lat.x", "min_lat.y")
#
# cbind(c, max_lat, min_lat)
# }, cl = cl)
#
# stopCluster(cl)
#
# centroids <- do.call("rbind",centroids)
# centroids$ID <- fc_list
#
# write.csv(centroids, "Data/centroids_study_areas.csv", row.names = FALSE)
centroids <- read.csv(here::here("Data/centroids_study_areas.csv"))
biov1 <- merge(biov1, centroids, by = "ID")
Following the approach proposed by Romain Bertrand
Just load it in
biov1 <- read.csv("Data/biov1_corrected_shifts.csv")
Identify direction of shift
biov1 <- biov1 %>%
mutate(shift_sign = ifelse(SHIFT>0,"pos","neg"),
shift_vel_sign = paste0(shift_sign,vel_sign),
shiftC_sign = ifelse(SHIFT_cor>0,"pos","neg"),
shiftC_vel_sign = paste0(shiftC_sign,vel_sign)
) %>%
filter(!is.na(velocity))
ggplot(biov1, aes(x=shift_vel_sign))+
geom_bar()+
coord_flip()+
facet_wrap(Type~Param)
ggplot(biov1, aes(x=shiftC_vel_sign))+
geom_bar()+
coord_flip()+
facet_wrap(Type~Param)
summary(biov1$velocity)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -5.8328 0.5415 1.9642 2.2647 2.6421 14.5492
summary(biov1$velocity[which(biov1$shift_vel_sign=="pospos")])
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.009396 0.541507 2.055833 2.462087 2.831932 14.549230
summary(biov1$velocity[which(biov1$shift_vel_sign=="negneg")])
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -3.79347 -2.30487 -1.78004 -1.60769 -0.58152 -0.02945
# calculate lags
# positive values are a lag (range shift lower smaller than expected or in opposite directions) and negative values are range shift larger than expected
biov1$lag <- biov1$velocity-biov1$SHIFT
position = which(biov1$vel_sign == "neg")
biov1$lag[position] <- biov1$lag[position] * -1 # Any negative velocity means the sign of lag has to shift.
# Lag 2
# When shift is in opposite sign of velocity, lag = abs(velocity)
biov1$lag2 = biov1$lag
position <- which(biov1$shift_vel_sign == "posneg" | biov1$shift_vel_sign == "negpos")
biov1$lag2[position] <- abs(biov1$velocity[position])
# Lag 3
# When shift is > velocity, lag = 0
biov1$lag3 = biov1$lag2
position <- which((biov1$shift_sign == "pos" & biov1$SHIFT > biov1$velocity) |
(biov1$shift_sign == "neg" & biov1$SHIFT < biov1$velocity))
biov1$lag3[position] <- 0
{
par(mfrow=c(2,2))
plot(lag~velocity, biov1)
plot(lag2~velocity, biov1)
plot(lag3~velocity, biov1)
}
ggplot(biov1, aes(x = Param, y=lag, fill = Param, color = Param))+
geom_boxplot(alpha = .5, outlier.shape = NA)+
scale_y_continuous(limits = quantile(biov1$lag2, c(0.1, 0.9), na.rm = T))+
facet_wrap(.~Type, scales = "free")
## Warning: Removed 8756 rows containing non-finite values (`stat_boxplot()`).
ggplot(biov1, aes(x = Param, y=lag2, fill = Param, color = Param))+
geom_boxplot(alpha = .5, outlier.shape = NA)+
scale_y_continuous(limits = quantile(biov1$lag2, c(0.1, 0.9), na.rm = T))+
facet_wrap(.~Type, scales = "free")
## Warning: Removed 6006 rows containing non-finite values (`stat_boxplot()`).
ggplot(biov1, aes(x = Param, y=log1p(lag3), fill = Param, color = Param))+
geom_boxplot(alpha = .5, outlier.shape = NA)+
# scale_y_continuous(limits = quantile(log1p(biov1$lag3), c(0.1, 0.9), na.rm = T))+
facet_wrap(.~Type, scales = "free")
# calculate lagCs
# positive values are a lagC (range shift lower smaller than expected or in opposite directions) and negative values are range shift larger than expected
biov1$lagC <- biov1$velocity-biov1$SHIFT_cor
position = which(biov1$vel_sign == "neg")
biov1$lagC[position] <- biov1$lagC[position] * -1 # Any negative velocity means the sign of lagC has to shift.
# Lag 2
# When shift is in opposite sign of velocity, lagC == velocity
biov1$lagC2 = biov1$lagC
position <- which(biov1$shiftC_vel_sign == "posneg" | biov1$shiftC_vel_sign == "negpos")
biov1$lagC2[position] <- abs(biov1$velocity[position])
# Lag 3
# When shift is > velocity, lagC = 0
biov1$lagC3 = biov1$lagC2
position <- which((biov1$shiftC_sign == "pos" & biov1$SHIFT_cor > biov1$velocity) |
(biov1$shiftC_sign == "neg" & biov1$SHIFT_cor < biov1$velocity))
biov1$lagC3[position] <- 0
{
par(mfrow=c(2,2))
plot(lagC~velocity, biov1)
plot(lagC2~velocity, biov1)
plot(lagC3~velocity, biov1)
}
ggplot(biov1, aes(x = Param, y=lagC, fill = Param, color = Param))+
geom_boxplot(alpha = .5, outlier.shape = NA)+
scale_y_continuous(limits = quantile(biov1$lagC2, c(0.1, 0.9), na.rm = T))+
facet_wrap(.~Type, scales = "free")
## Warning: Removed 6070 rows containing non-finite values (`stat_boxplot()`).
ggplot(biov1, aes(x = Param, y=lagC2, fill = Param, color = Param))+
geom_boxplot(alpha = .5, outlier.shape = NA)+
scale_y_continuous(limits = quantile(biov1$lagC2, c(0.1, 0.9), na.rm = T))+
facet_wrap(.~Type, scales = "free")
## Warning: Removed 6012 rows containing non-finite values (`stat_boxplot()`).
ggplot(biov1, aes(x = Param, y=log1p(lagC3), fill = Param, color = Param))+
geom_boxplot(alpha = .5, outlier.shape = NA)+
# scale_y_continuous(limits = quantile(log1p(biov1$lagC3), c(0.1, 0.9), na.rm = T))+
facet_wrap(.~Type, scales = "free")
plot(SHIFT_cor~SHIFT,biov1)
abline(a=0,b=1,col=2)
lm0=lm(SHIFT_cor~SHIFT,biov1)
summary(lm0)
##
## Call:
## lm(formula = SHIFT_cor ~ SHIFT, data = biov1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -19.6130 -0.1863 0.0470 0.1971 10.9133
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.1815049 0.0028742 -63.15 <2e-16 ***
## SHIFT 0.1556612 0.0005166 301.33 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4872 on 30055 degrees of freedom
## Multiple R-squared: 0.7513, Adjusted R-squared: 0.7513
## F-statistic: 9.08e+04 on 1 and 30055 DF, p-value: < 2.2e-16
hist(biov1$SHIFT,col=rgb(1,0,0,0.5))
hist(biov1$SHIFT_cor,col=rgb(0,0,1,0.5),add=T)
#so by using residuals we decrease the variation in the raw range shift obs, it's like we work on a more homogenous variables as all the variations due to methods variation has been substracted to the shift.
x1=tapply(biov1$SHIFT,biov1$Group,mean)
x2=tapply(biov1$SHIFT_cor,biov1$Group,mean)#lower mean value than true obs
plot(x1~x2, xlab="Mean shift per Group (raw)", ylab="Mean shift per Group (corrected)")
lm0=lm(x1~x2,biov1)
summary(lm0) #it changes many things
##
## Call:
## lm(formula = x1 ~ x2, data = biov1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.5194 -1.0516 -0.2973 0.4605 3.3886
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.1183 0.3504 6.045 2.24e-05 ***
## x2 5.6218 2.7006 2.082 0.0549 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.418 on 15 degrees of freedom
## Multiple R-squared: 0.2241, Adjusted R-squared: 0.1724
## F-statistic: 4.333 on 1 and 15 DF, p-value: 0.05491
x1=tapply(biov1$SHIFT,biov1$Group,var)
x2=tapply(biov1$SHIFT_cor,biov1$Group,var)#lower variance value than in true obs
plot(x1~x2, xlab="Variance shift per Group (raw)", ylab="Variance shift per Group (corrected)")
lm0=lm(x1~x2,biov1)
summary(lm0) #but high positive correlation meaning that relative variation is conserved
##
## Call:
## lm(formula = x1 ~ x2, data = biov1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -21.801 -19.391 -9.706 -1.438 71.182
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.755 8.042 2.829 0.0127 *
## x2 17.645 2.262 7.801 1.17e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 29.41 on 15 degrees of freedom
## Multiple R-squared: 0.8022, Adjusted R-squared: 0.7891
## F-statistic: 60.85 on 1 and 15 DF, p-value: 1.173e-06
x1=tapply(biov1$SHIFT,biov1$sp_name_std_v1,mean)
x2=tapply(biov1$SHIFT_cor,biov1$sp_name_std_v1,mean)#lower mean value than in true obs
plot(x1~x2, xlab="Mean shift per species (raw)", ylab="Mean shift per species (corrected)")
lm0=lm(x1~x2,biov1)
summary(lm0) #but at the species level, we observe a high positive relationship so the corrected shifts did not change the pattern of range shift=> good
##
## Call:
## lm(formula = x1 ~ x2, data = biov1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -32.809 -1.107 -0.470 0.569 129.850
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.41199 0.02744 51.47 <2e-16 ***
## x2 5.04851 0.03134 161.08 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.02 on 12118 degrees of freedom
## Multiple R-squared: 0.6817, Adjusted R-squared: 0.6816
## F-statistic: 2.595e+04 on 1 and 12118 DF, p-value: < 2.2e-16
#looking at the relationship with climate velocity
#WARNING:here you need to adapt the code: if you look at the elevational range shift so you have to use EleVeloT in order to use the corresponding climate velocity
par(mfrow=c(1,2))
plot(SHIFT~velocity, xlab="SHIFT (raw)", ylab="Velocity",biov1)
lm0=lm(SHIFT~velocity+I(velocity^2),biov1)
summary(lm0)
##
## Call:
## lm(formula = SHIFT ~ velocity + I(velocity^2), data = biov1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -120.794 -1.642 -0.792 1.232 145.242
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.849929 0.049121 17.303 < 2e-16 ***
## velocity 0.195422 0.028720 6.804 1.03e-11 ***
## I(velocity^2) -0.010881 0.002856 -3.810 0.000139 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 5.434 on 30054 degrees of freedom
## Multiple R-squared: 0.002526, Adjusted R-squared: 0.00246
## F-statistic: 38.06 on 2 and 30054 DF, p-value: < 2.2e-16
x1=seq(min(biov1$velocity,na.rm = T),max(biov1$velocity,na.rm = T),le=100)
p1=predict(lm0,newdata=data.frame(velocity=x1),type="response")
points(p1~x1,type="l",col=2,lwd=2)
plot(SHIFT_cor~velocity, xlab="SHIFT (corrected)", ylab="Velocity",biov1)#the pattern is changing a lot
lm0=lm(SHIFT_cor~velocity+I(velocity^2),biov1)
summary(lm0)
##
## Call:
## lm(formula = SHIFT_cor ~ velocity + I(velocity^2), data = biov1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -18.1763 -0.3541 -0.0312 0.2835 18.5895
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.0091560 0.0088312 -1.037 0.29985
## velocity 0.0115201 0.0051635 2.231 0.02568 *
## I(velocity^2) -0.0014568 0.0005135 -2.837 0.00456 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9769 on 30054 degrees of freedom
## Multiple R-squared: 0.0002878, Adjusted R-squared: 0.0002213
## F-statistic: 4.326 on 2 and 30054 DF, p-value: 0.01322
x1=seq(min(biov1$velocity,na.rm = T),max(biov1$velocity,na.rm = T),le=100)
p1=predict(lm0,newdata=data.frame(velocity=x1),type="response")
points(p1~x1,type="l",col=2,lwd=2)
#we observed a relationship for the raw range shift (an unimodal relationship with higher shifts when absolute climate velocity increase).
#for the range shift corrected by method, no relationship with climatic velocity is observed
#looking at the relationship with lag
par(mfrow=c(1,2))
plot(lagC~velocity,biov1,main="corrected lag estimate")
# Lag C2
# When shift is in opposite sign of velocity, lag == velocity
plot(lagC2~velocity,biov1,main="corrected lag estimate with correction for special case")
#change more things than the above lag metrics, still the highly change are observed at extreme negative values
plot(lag2~lagC2,biov1)
lm0=lm(lag~lagC,biov1)
summary(lm0) #we observe a high positive relationship so the corrected lag did not change the pattern of the observed lag=> good
##
## Call:
## lm(formula = lag ~ lagC, data = biov1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -140.499 -1.276 0.724 1.869 96.858
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.056745 0.034810 -59.09 <2e-16 ***
## lagC 1.414888 0.009807 144.27 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.509 on 30055 degrees of freedom
## Multiple R-squared: 0.4092, Adjusted R-squared: 0.4091
## F-statistic: 2.081e+04 on 1 and 30055 DF, p-value: < 2.2e-16
par(mfrow=c(1,2))
plot(lag2~velocity,biov1)
lm0=lm(lag2~velocity+I(velocity^2),biov1)
summary(lm0) #significant; R2=46%
##
## Call:
## lm(formula = lag2 ~ velocity + I(velocity^2), data = biov1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -144.671 -0.360 1.098 2.100 6.457
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.874182 0.037989 -23.01 <2e-16 ***
## velocity 0.273585 0.022211 12.32 <2e-16 ***
## I(velocity^2) 0.054244 0.002209 24.56 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.202 on 30054 degrees of freedom
## Multiple R-squared: 0.191, Adjusted R-squared: 0.1909
## F-statistic: 3548 on 2 and 30054 DF, p-value: < 2.2e-16
x1=seq(min(biov1$velocity, na.rm = T),max(biov1$velocity, na.rm = T),le=100)
p1=predict(lm0,newdata=data.frame(velocity=x1),type="response")
points(p1~x1,type="l",col=2,lwd=2) #bimodal relationship
plot(lagC2~velocity,biov1)
lm0=lm(lagC2~velocity+I(velocity^2),biov1)
summary(lm0)#significant; R2=86
##
## Call:
## lm(formula = lagC2 ~ velocity + I(velocity^2), data = biov1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -18.4839 -0.2248 0.0613 0.3711 7.2808
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.2523552 0.0075364 33.48 <2e-16 ***
## velocity 0.6480167 0.0044064 147.06 <2e-16 ***
## I(velocity^2) 0.0305757 0.0004382 69.78 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.8337 on 30054 degrees of freedom
## Multiple R-squared: 0.891, Adjusted R-squared: 0.8909
## F-statistic: 1.228e+05 on 2 and 30054 DF, p-value: < 2.2e-16
x1=seq(min(biov1$velocity, na.rm = T),max(biov1$velocity, na.rm = T),le=100)
p1=predict(lm0,newdata=data.frame(velocity=x1),type="response")
points(p1~x1,type="l",col=2,lwd=2) #bimodal relationship
ggplot(biov1, aes(x = velocity, y = SHIFT))+
geom_point(aes(color = lag2), alpha = .1)+
scale_color_viridis_c()+
geom_smooth(method = "lm")+
theme_classic()+
facet_wrap(Param~Type+ECO, scales = "free")
## `geom_smooth()` using formula = 'y ~ x'
ggplot(biov1 %>%
filter(shift_vel_sign == "pospos" | shift_vel_sign == "negneg"),
aes(x = abs(velocity), y = abs(SHIFT)))+
geom_point(aes(color = lag2), alpha = .1)+
scale_color_viridis_c()+
geom_smooth(method = "lm")+
theme_classic()+
facet_wrap(Param~Type+ECO, scales = "free")
## `geom_smooth()` using formula = 'y ~ x'
ggplot(biov1, aes(x = velocity, y = SHIFT_cor))+
geom_point(aes(color = lag2), alpha = .1)+
scale_color_viridis_c()+
geom_smooth(method = "lm")+
theme_classic()+
facet_wrap(Param~Type+ECO, scales = "free")
## `geom_smooth()` using formula = 'y ~ x'
ggplot(biov1 %>%
filter(shift_vel_sign == "pospos" | shift_vel_sign == "negneg"),
aes(x = abs(velocity), y = abs(SHIFT_cor)))+
geom_point(aes(color = lag2), alpha = .1)+
scale_color_viridis_c()+
geom_smooth(method = "lm")+
theme_classic()+
facet_wrap(Param~Type+ECO, scales = "free")
## `geom_smooth()` using formula = 'y ~ x'
# Genetic diversity
mt <- fread(here::here('Data/mtdna.csv'))
mt$spp <- Clean_Names(mt$spp)
ms <- fread(here::here('Data/msat.csv'))
ms$spp <- Clean_Names(ms$spp)
gen_d <- fread(here::here('Data/Deposited_data_genetic_diversity_dekort2021.csv'))
gen_d$spp <- gsub("_"," ",gen_d$Species)
gen_d$spp <- Clean_Names(gen_d$spp)
gen_lf <- read_excel(here::here("Data/MacroPopGen_Database_final-v0.2.xlsx"), sheet = "MacroPopGen_Database")
gen_lf$spp <- gsub("_"," ",gen_lf$G_s)
gen_lf$spp <- Clean_Names(gen_lf$spp)
# Phenotypic rates
pheno <- fread(here::here('Data/PROCEEDv6_RatesDB.csv'))
pheno$spp <- Clean_Names(pheno$sp_ncbi)
length(unique(gen_d$spp))
length(unique(gen_lf$spp))
length(unique(mt$spp))
length(unique(ms$spp))
length(unique(pheno$spp))
# Species to find
mycols <- c("species","scientificName","kingdom","phylum","class","order","family","db_code")
splist <- data.frame(matrix(ncol = length(mycols), nrow = length(unique(c(gen_sps,phen_sps)))))
names(splist) <- mycols
splist$reported_name_fixed <- unique(c(gen_sps,phen_sps))
tofind_ <- splist[which(is.na(splist$scientificName)),]
tofind_ <- unique(tofind_$reported_name_fixed)
tofind <- data.frame(matrix(nrow = length(tofind_), ncol = 8))
names(tofind) = c("scientificName", "kingdom", "phylum", "class", "order", "family", "db", "db_code")
tofind <- data.frame(species = tofind_, tofind)
tofind <- tofind %>%
mutate(across(everything(), as.character))
# retrieve sp names
sp_names_found <- Find_Sci_Names(spnames = tofind$species)
# Feed
tofind <- tofind %>%
rows_patch(sp_names_found,
by = "species")
gc()
## Add found species names to the splist
all(tofind$species %in% splist$reported_name_fixed)
if(any(is.na(tofind$scientificName))){
found <- tofind[-which(is.na(tofind$scientificName)),]
} else {
found = tofind
}
for(i in 1:length(found$species)){
tofill <- unique(which(splist$reported_name_fixed == found$species[i]))
splist$scientificName[tofill] <- found$scientificName[i]
splist$kingdom[tofill] <- found$kingdom[i]
splist$phylum[tofill] <-found$phylum[i]
splist$class[tofill] <- found$class[i]
splist$order[tofill] <- found$order[i]
splist$family[tofill] <- found$family[i]
splist$db[tofill] <- found$db[i]
splist$db_code[tofill] <- found$db_code[i]
}
splist$spp = splist$reported_name_fixed
# Keep original names for those species we could not find a name at GBIF
pos <- which(is.na(splist$scientificName))
splist$scientificName[pos] <- splist$spp[pos]
any(!mt$spp %in% gen_sps)
any(!mt$spp %in% splist$spp)
any(!ms$spp %in% gen_sps)
any(!ms$spp %in% splist$spp)
any(!gen_d$spp %in% gen_sps)
any(!gen_d$spp %in% splist$spp)
any(!gen_lf$spp %in% gen_sps)
any(!gen_lf$spp %in% splist$spp)
for(i in 1:nrow(mt)){
mt$spp_new[i] <- splist$scientificName[which(splist$spp == mt$spp[i])]
}
for(i in 1:nrow(ms)){
ms$spp_new[i] <- splist$scientificName[which(splist$spp == ms$spp[i])]
}
for(i in 1:nrow(gen_d)){
gen_d$spp_new[i] <- splist$scientificName[which(splist$spp == gen_d$spp[i])]
}
for(i in 1:nrow(gen_lf)){
gen_lf$spp_new[i] <- splist$scientificName[which(splist$spp == gen_lf$spp[i])]
}
for(i in 1:nrow(gen_lf)){
pheno$spp_new[i] <- splist$scientificName[which(splist$spp == pheno$spp[i])]
}
write.csv(mt,
here::here('Data/mtdna_harmo.csv'), row.names = FALSE)
write.csv(ms,
here::here('Data/msat_harmo.csv'), row.names = FALSE)
write.csv(gen_d,
here::here('Data/Deposited_data_genetic_diversity_dekort2021_harmo.csv'), row.names = FALSE)
write.csv(gen_lf,
here::here('Data/MacroPopGen_Database_final_areas_Lawrence_Fraser_2020_harmo.csv'), row.names = FALSE)
write.csv(pheno,
here::here('Data/PROCEEDv6_RatesDB_harmo.csv'), row.names = FALSE)
mt <- fread(here::here('Data/mtdna_harmo.csv'))
ms <- fread(here::here('Data/msat_harmo.csv'))
gen_d <- fread(here::here('Data/Deposited_data_genetic_diversity_dekort2021_harmo.csv'))
gen_d$spp <- gsub("_"," ",gen_d$Species)
gen_lf <- fread(here::here('Data/MacroPopGen_Database_final_areas_Lawrence_Fraser_2020_harmo.csv'))
gen_lf$spp <- gsub("_"," ",gen_lf$G_s)
pheno <- fread(here::here('Data/PROCEEDv6_RatesDB_harmo.csv'))
pheno$spp <- gsub("_"," ",pheno$sp_ncbi)
# Gen
gen_sps <- unique(c(mt$spp, ms$spp, gen_d$spp, gen_lf$spp))
bsi_v1 <- intersect(biov1$spp, gen_sps)
# Break down
tot1 <- data.frame(N_species_with_gen_data = length(gen_sps),
v1_matches = length(bsi_v1))
tot1
## N_species_with_gen_data v1_matches
## 1 1894 287
# phenotypic
phen_sps <- unique(pheno$spp)
bsj_v1 <- intersect(biov1$spp, phen_sps)
# Break down
tot2 <- data.frame(N_species_with_pheno_data = length(phen_sps),
v1_matches = length(bsj_v1))
tot2
## N_species_with_pheno_data v1_matches
## 1 428 260
# Malin
ggplot(mt, aes(He)) +
geom_histogram() +
ggtitle("Malin's Mitochondrial")
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## Warning: Removed 48 rows containing non-finite values (`stat_bin()`).
ggplot(ms, aes(He)) +
geom_histogram() +
ggtitle("Malin's Microsatellite")
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
# De Kort
ggplot(gen_d, aes(GDp)) +
geom_histogram() +
ggtitle("De Kort AFLP & Microsatellite")
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
# Why De Kort data has the two peaks? It is due to the marker type
ggplot(gen_d, aes(GDp)) +
geom_histogram()+
facet_grid(.~Marker)
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
# But if we look at the normalized GDp, the problem is solved. Therefore, we should be using the normalized version of GDp from DeKort
ggplot(gen_d, aes(GDp_norm)) +
geom_histogram()
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
# Lawrence & Fraser
ggplot(gen_lf, aes(as.numeric(He))) +
geom_histogram()+
ggtitle("Lawrence & Fraser Microsatellite")
## Warning in FUN(X[[i]], ...): NAs introduced by coercion
## Warning in FUN(X[[i]], ...): NAs introduced by coercion
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## Warning: Removed 1283 rows containing non-finite values (`stat_bin()`).
# Malin
mt_sub <- mt %>%
mutate(long = as.numeric(lon),
lat = as.numeric(lat),
N = as.numeric(n), # Sample size
Marker = "Mitochondrial") %>%
dplyr::filter(!is.na(He) | !He==0) %>%
dplyr::select(spp, He, long, lat, Marker)
mt_sub$Source = "Malin"
mt_sub$Het_type = "He"
ms_sub <- ms %>%
mutate(long = as.numeric(lon),
lat = as.numeric(lat),
N = as.numeric(n), # Sample size
Marker = "Microsatellite") %>%
dplyr::filter(!is.na(He) | !He==0) %>%
dplyr::select(spp, He, long, lat, Marker)
ms_sub$Source = "Malin"
ms_sub$Het_type = "He"
###################
# De Kort
# Jonathan R.:
# co-dominant refers to the fact that you can distinguish homozigous and heterozygous individuals, such as microsatellites
# dominant means that you either have the detection of the marker or not. but you don’t know if the individual is heterozygous or not, such as AFLP
# restriction enzymes are proteins that you use to cut the genomes into fragments. there are used for AFLP and allozymes, so not sure
# De Kort code for Marker:
# CD = Co-dominant
# D = Dominant
# ENZ = Enzime
gen_d_sub <- gen_d %>%
mutate(long = as.numeric(LONGITUDE),
lat = as.numeric(LATITUDE),
spp = spp_new,
He = as.numeric(GDp),
N = as.numeric(SampleSize), # Sample size
Marker = case_when(
Marker=="CD" ~ "Microsatellite",
Marker=="D" ~ "AFLP",
Marker=="ENZ" ~ "AFLP",
TRUE ~ as.character(Marker))) %>%
dplyr::filter(!is.na(He) | !He==0) %>%
dplyr::filter(!(Marker == "AFLP" & He > .5)) %>% # AFLPs > 0.5 are errors
dplyr::select(spp, He, long, lat, Marker)
gen_d_sub$Source = "De Kort"
gen_d_sub$Het_type = "He"
# How DeKort calculated the normalized version of GDp
# gen_d_sub$He_harm <- NA
#
# m <- unique(gen_d_sub$Marker)
# for(i in 1:length(m)){
# pos <- which(gen_d_sub$Marker == m[i])
# gen_d_sub$He_harm[pos] <- scale(gen_d_sub$He[pos])
# }
# gen_d_sub$He_harm <- (gen_d_sub$He_harm - min(gen_d_sub$He_harm))/(max(gen_d_sub$He_harm)-min(gen_d_sub$He_harm))
###################
# Lawrence & Fraser
pos <- which(is.na(gen_lf$He))
gen_lf$He_new <- gen_lf$He
gen_lf$He_new[pos] <- gen_lf$Ho[pos]
gen_lf$Het_type = "He"
gen_lf$Het_type[pos] = "Ho"
gen_lf_sub <- gen_lf %>%
mutate(long = as.numeric(Long),
lat = as.numeric(Lat),
He = as.numeric(He_new),
spp = spp_new,
N = as.numeric(n), # Sample size
Marker = "Microsatellite") %>%
dplyr::filter(!is.na(He) | !He==0) %>%
dplyr::select(spp, He, long, lat, Marker, Het_type)
## Warning: There were 2 warnings in `mutate()`.
## The first warning was:
## ℹ In argument: `He = as.numeric(He_new)`.
## Caused by warning:
## ! NAs introduced by coercion
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 1 remaining warning.
gen_lf_sub$Source = "Lawrence & Fraser"
###################
# Group all
gen_div <- rbind(mt_sub %>% dplyr::select(names(mt_sub)),
ms_sub %>% dplyr::select(names(mt_sub)),
gen_d_sub %>% dplyr::select(names(mt_sub)),
gen_lf_sub %>% dplyr::select(names(mt_sub)))
gen_div <- na.omit(gen_div)
gen_div$spp <- as.factor(gen_div$spp)
gen_div$Marker <- as.factor(gen_div$Marker)
gen_div$Source <- as.factor(gen_div$Source)
gen_div$He_harm <- NA
gen_div$He_harm2 <- NA
gen_div$He_harm3 <- NA
m <- unique(gen_div$Marker)
for(i in 1:length(m)){
# Select marker type
pos <- which(gen_div$Marker == m[i])
# Apply DeKort transformation
gen_div$He_harm[pos] <- deKort_trans(gen_div$He[pos])
# Apply logit transformation
gen_div$He_harm2[pos] <- logit_trans(gen_div$He[pos])
# Centered in zero
gen_div$He_harm3[pos] <- scale(gen_div$He_harm2[pos])
}
# Filter the species in Bioshifts
gen_div <- gen_div %>% dplyr::filter(spp %in% unique(biov1$spp))
# remove outliers
gen_metrics <- c("He_harm","He_harm2","He_harm3")
for (i in 1:length(gen_metrics)){
my_gem <- as.numeric(data.frame(gen_div %>% select(gen_metrics[i]))[,1])
quartiles <- quantile(my_gem, probs=c(.05, .95), na.rm = FALSE)
IQR <- IQR(my_gem)
Lower <- quartiles[1] - 1.5*IQR
Upper <- quartiles[2] + 1.5*IQR
gen_div <- subset(gen_div, my_gem > Lower & my_gem < Upper)
}
We averaged genetic diversity by location (same species, in the same XY coordinates)
# How many species with >1 gen div measurements at the same site?
# Check how many obs per site and marker type exist
N_obs <- gen_div %>%
group_by(spp, long, lat) %>%
tally() %>%
dplyr::filter(n>1)
DT::datatable(N_obs)
# Is there any species with >1 marker type in the same location?
gen_div %>%
group_by(spp, long, lat, Marker) %>%
dplyr::summarise(N = length(unique(spp))) %>%
dplyr::filter(N>1)
## `summarise()` has grouped output by 'spp', 'long', 'lat'. You can override
## using the `.groups` argument.
## # A tibble: 0 × 5
## # Groups: spp, long, lat [0]
## # ℹ 5 variables: spp <fct>, long <dbl>, lat <dbl>, Marker <fct>, N <int>
# No, there are not
# Is there any species with >1 Source type in the same location?
gen_div %>%
group_by(spp, long, lat, Source) %>%
dplyr::summarise(N = length(unique(spp))) %>%
dplyr::filter(N>1)
## `summarise()` has grouped output by 'spp', 'long', 'lat'. You can override
## using the `.groups` argument.
## # A tibble: 0 × 5
## # Groups: spp, long, lat [0]
## # ℹ 5 variables: spp <fct>, long <dbl>, lat <dbl>, Source <fct>, N <int>
# No, there are not
# Avg by site and species
gen_div <- gen_div %>%
group_by(spp, long, lat) %>%
dplyr::summarise(He = median(He),
He_harm = median(He_harm),
He_harm2 = median(He_harm2),
He_harm3 = median(He_harm3),
Source = paste(Source, sep = ", "),
Marker = paste(Marker, sep = ", "))
## Warning: Returning more (or less) than 1 row per `summarise()` group was deprecated in
## dplyr 1.1.0.
## ℹ Please use `reframe()` instead.
## ℹ When switching from `summarise()` to `reframe()`, remember that `reframe()`
## always returns an ungrouped data frame and adjust accordingly.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## `summarise()` has grouped output by 'spp', 'long', 'lat'. You can override
## using the `.groups` argument.
unique(gen_div$Source)
## [1] "De Kort" "Lawrence & Fraser" "Malin"
unique(gen_div$Marker)
## [1] "AFLP" "Microsatellite" "Mitochondrial"
gen_div <- merge(gen_div,
splist[,-7],
by.x="spp",
by.y="scientificName",
all.x = T)
toplot_ <- gen_div %>%
group_by(Group) %>%
dplyr::summarise(Obs = length(spp),
Species = length(unique(spp))) %>%
gather(var, N, -c(Group))
ggplot(toplot_, aes(x = Group, y = N))+
ggtitle("N species by Group")+
geom_bar(stat="identity")+
geom_text(aes(y=N, label=N,
hjust = ifelse(N<max(N),-.1,1.1)), vjust=0.2, size=3,
position = position_dodge(0.9))+
theme_classic()+
coord_flip()+
facet_wrap(.~var, scales = "free", ncol=1)
# correlations
ggpairs(gen_div,
columns = c("He","He_harm","He_harm2","He_harm3"))
ggpairs(gen_div,
columns = c("He","He_harm","He_harm2","He_harm3"),
aes(color = Marker,
alpha = 0.5))
ggpairs(gen_div,
columns = c("He","He_harm","He_harm2","He_harm3"),
aes(color = Source,
alpha = 0.5))
to_plot <- gen_div %>%
select(Group,He_harm,He_harm2,He_harm3) %>%
filter(Group %in% c("Mammal","Fish","Plant", "Bird")) %>%
gather(metric, value, -Group)
ggplot(to_plot, aes(value, fill = Group, color = Group))+
geom_density(alpha = .3)+
scale_color_viridis_d()+
scale_fill_viridis_d()+
facet_wrap(.~metric, scales = "free")
toplot <- gen_div %>%
mutate(Hemisphere = ifelse(lat<0, "South", "North"),
lat_abs = abs(lat),
He_harm = scale(He_harm),
He_harm3 = scale(He_harm3))
ggplot(toplot, aes(x = lat_abs, y = He_harm, color = Group))+
geom_point(alpha = .1)+
stat_smooth(method = "lm")+
facet_wrap(.~Group)
## `geom_smooth()` using formula = 'y ~ x'
mod1 <- lm(He_harm~lat*Group, toplot)
emtrends(mod1, ~ lat*Group, var = "lat")
## lat Group lat.trend SE df lower.CL upper.CL
## 46.9 Amphibian -0.06867 0.029720 7833 -0.12693 -0.01041
## 46.9 Bird 0.00394 0.001394 7833 0.00121 0.00668
## 46.9 Fish 0.01032 0.000837 7833 0.00868 0.01196
## 46.9 Mammal -0.00464 0.004656 7833 -0.01377 0.00449
## 46.9 Molluscs 0.07562 0.186554 7833 -0.29008 0.44131
## 46.9 Plant -0.01902 0.001150 7833 -0.02128 -0.01677
## 46.9 Reptile 0.08038 0.023489 7833 0.03433 0.12642
##
## Confidence level used: 0.95
mundi <- terra::vect(rnaturalearth::ne_coastline(scale = 110, returnclass = "sp"))
## Warning: The `returnclass` argument of `ne_download()` sp as of rnaturalearth 1.0.0.
## ℹ Please use `sf` objects with {rnaturalearth}, support for Spatial objects
## (sp) will be removed in a future release of the package.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
# get vect gen div
vect_div <- terra::vect(gen_div,geom=c("long","lat"),crs=crs(mundi))
# get vect centroids bioshifts
vect_biov1 <- terra::vect(biov1 %>%
filter(spp %in% gen_div$spp),
geom=c("centroid.x","centroid.y"),crs=crs(mundi))
# get raster bioshifts shp files study areas
get_raster_bioshifts = "NO"
if(get_raster_bioshifts=="YES"){
my_ext = terra::ext(mundi)
my_crs = crs(mundi)
rast_biov1 <- terra::rast(my_ext, crs = my_crs, res = 0.5)
values(rast_biov1) <- 0
fgdb <- "C:/Users/brunn/NextCloud/Bioshifts/Study_Areas.gdb"
fc_list <- ogrListLayers(fgdb)
fc_list <- fc_list[which(fc_list %in% unique(vect_biov1$ID))]
for(i in 1:length(fc_list)){ cat("\r",i,"from",length(fc_list))
tmp = terra::vect(sf::st_read(fgdb, layer=fc_list[i]))
tmp = terra::cells(rast_biov1,tmp)
tmp_cell = tmp[,2]
tmp_vals = rast_biov1[tmp_cell][,1]
rast_biov1[tmp_cell] = tmp_vals+1
}
names(rast_biov1) <- "SA"
rast_biov1[rast_biov1==0] <- NA
writeRaster(rast_biov1, "Data/raster_bioshifts_SA.tif")
} else {
rast_biov1 <- terra::rast("Data/raster_bioshifts_SA.tif")
}
values_biov1 <- na.omit(values(rast_biov1))
ggplot()+
ggtitle("Genetic diversity data")+
geom_spatvector(data=mundi)+
geom_spatvector(data=vect_div,aes(color = Marker))+
theme_blank()
ggplot()+
ggtitle("Genetic diversity data")+
geom_spatvector(data=mundi)+
geom_spatvector(data=vect_div,aes(color = Source))+
theme_blank()
ggplot()+
ggtitle("Bioshifts data \n(Centroid of study areas)")+
geom_spatvector(data=mundi)+
geom_spatvector(data=vect_biov1, aes(color = ECO))+
theme_blank()
ggplot()+
ggtitle("Bioshifts data \n(Overlap study areas)")+
geom_spatraster(data=rast_biov1)+
scale_fill_whitebox_b(
palette = "muted",
na.value = "white",
breaks = seq(min(values_biov1), max(values_biov1), 1))+
geom_spatvector(data=vect_biov1, aes(color = ECO))+
geom_spatvector(data=mundi)+
theme_blank()
pheno_rate <- pheno %>%
mutate(long = as.numeric(sample1.longitude),
lat = as.numeric(sample1.latitude),
trait1 = as.numeric(mean1),
trait2 = as.numeric(mean2),
spp = spp_new,
Trait = trait_type,
Years = as.numeric(years),
Gen = as.numeric(generations)) %>%
dplyr::select(spp, long, lat, Trait, trait1, trait2, Years, Gen) %>%
dplyr::filter(!is.na(trait1) & !is.na(trait2))
pheno_rate$Source = "PROCEED"
for(i in 1:nrow(pheno_rate)){
sppi <- pheno_rate$spp[i]
traiti <- pheno_rate$Trait[i]
pos_spp_trait_i <- which(pheno_rate$spp == sppi & pheno_rate$Trait == traiti)
pheno_rate$PR1[i] <- (pheno_rate$trait2[i] - pheno_rate$trait1[i]) / mean(c(pheno_rate$trait2[pos_spp_trait_i], pheno_rate$trait1[pos_spp_trait_i])) * pheno_rate$Years[i]
pheno_rate$PR2[i] <- (pheno_rate$trait2[i] - pheno_rate$trait1[i]) / sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i], pheno_rate$trait1[pos_spp_trait_i])))^2 * pheno_rate$Years[i]
}
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
## Warning in sqrt(mean(c(pheno_rate$trait2[pos_spp_trait_i],
## pheno_rate$trait1[pos_spp_trait_i]))): NaNs produced
pheno_rate$PR_harm <- NA
m <- unique(pheno_rate$Trait)
for(i in 1:length(m)){
pos <- which(pheno_rate$Trait == m[i])
# DeKort
p <- scale(pheno_rate$PR1[pos])
p <- (p - min(p, na.rm = T))/(max(p, na.rm = T)-min(p, na.rm = T))
pheno_rate$PR_harm[pos] <- p
}
# gen_div %>%
# group_by(spp) %>%
# tally() %>%
# summary()
# Look into phenotic rates
ggplot(pheno_rate, aes(PR1)) +
geom_histogram() +
facet_wrap(Trait~., scales = "free") +
ggtitle("Raw data")
## Warning: Removed 1 rows containing non-finite values (`stat_bin()`).
ggplot(pheno_rate, aes(PR_harm)) +
geom_histogram() +
facet_wrap(Trait~., scales = "free") +
ggtitle("Harmonization following De Kort")
## Warning: Removed 1 rows containing non-finite values (`stat_bin()`).
We tried several methods for merging genetic data and bioshifts.
Method 1: Normal merge Uses the merge function to aggregate bioshifts and gen div data. This creates multiple false duplicates.
Method 2: Median He per species. This ignores spatial proximity of He values to range shift values.
Method 3: Median He per Param & species Use lower latitude TE and higher latitude as LE (North Hemisphere). Although this collects He values based on latitude, the proximity of He to the study areas where range shifts were observed is ignored.
Method 4: Distance based approach. Calculates distance weighted mean of He to the centroid of the study area.
Method 5: Location based approach. Get the closest He to the centroid (for centroid shift), to the max latitude (for LE shift, at the North hemisphere), and to the min latitude (for TE shift, at the North hemisphere) of the study area.
We decided that Method 5 is the best one.
This is not the correct way to merge because it creates multiple false duplicates.
This ignores spatial proximity of He values to range shift values.
gen_div_avg <- gen_div %>%
group_by(spp) %>%
dplyr::summarise(N = length(spp),
He_harm = median(He_harm),
He_harm2 = median(He_harm2),
He_harm3 = median(He_harm3))
#merge all data frames in list
gen_data_v1_avg <- append(list(gen_div_avg), list(biov1)) %>% reduce(full_join, by='spp')
gen_data_v1_avg <- gen_data_v1_avg %>% dplyr::filter(!is.na(He_harm), !is.na(SHIFT))
# N species
length(unique(gen_data_v1_avg$spp))
## [1] 304
write.csv(gen_data_v1_avg,"Data/gen_data_final_m2.csv",row.names = FALSE)
Use lower latitude TE and higher latitude as LE (North
Hemisphere).
Although this collects He values based on latitude, the proximity of He
to the study areas where range shifts were observed is ignored.
Calculate distance weighted mean He for each species in the bioshifts
database.
This method collects He values more close to the centroid of the study
area.
Get the closest He to the centroid (for O), to the max latitude (for LE, at the North hemisphere), and to the min latitude (for TE, at the North hemisphere) of the study area.
# sp list
m <- unique(biov1$spp[which(biov1$spp %in% gen_div$spp)])
gen_data_v1_dist2 <- data.frame()
# sp2go="Alces alces"
# sp2go="Sylvia atricapilla"
for(i in 1:length(m)){ cat(i, "from", length(m), "\r")
sp2go <- m[i]
# subset genetic data
sub_gen <- subset(gen_div[,1:9], spp == sp2go)
sub_gen_v <- terra::vect(sub_gen,
geom=c("long", "lat"),
crs = "+proj=longlat +datum=WGS84 +no_defs")
# subset bioshifts data
sub_bio <- subset(biov1, spp == sp2go)
# study areas for species i
areas <- unique(sub_bio$ID)
gen_data_v1_dist2_tmp <- data.frame()
# get genetic data for species i in each study area
for(j in 1:nrow(sub_bio)){
# bioshifts in study area j and range pos r
sub_bio_area_j_tmp <- sub_bio[j,]
range_pos2go <- sub_bio_area_j_tmp$Param
coord_names <- NA
coord_names <- if(range_pos2go == "O"){ c("centroid.x", "centroid.y") } else {coord_names}
coord_names <- if(range_pos2go == "LE" & sub_bio_area_j_tmp$centroid.y > 0){ c("max_lat.x", "max_lat.y") } else {coord_names}
coord_names <- if(range_pos2go == "LE" & sub_bio_area_j_tmp$centroid.y < 0){ c("min_lat.x", "min_lat.y") } else {coord_names}
coord_names <- if(range_pos2go == "TE" & sub_bio_area_j_tmp$centroid.y < 0){ c("max_lat.x", "max_lat.y") } else {coord_names}
coord_names <- if(range_pos2go == "TE" & sub_bio_area_j_tmp$centroid.y > 0){ c("min_lat.x", "min_lat.y") } else {coord_names}
sub_bio_area_j_v = data.frame(sub_bio_area_j_tmp)[,coord_names]
names(sub_bio_area_j_v) <- c("long", "lat")
sub_bio_area_j_v <- terra::vect(unique(sub_bio_area_j_v[,c("long", "lat")]),
geom=c("long", "lat"),
crs = "+proj=longlat +datum=WGS84 +no_defs")
# plot(terra::vect(rbind(sub_gen[,c("long", "lat")], sub_bio_area_j_tmp[,c("long", "lat")]),
# geom=c("long", "lat"),
# crs = "+proj=longlat +datum=WGS84 +no_defs"))
# plot(sub_bio_area_j_v, col = "red", add=T)
dists <- terra::distance(sub_bio_area_j_v, sub_gen_v)[1,]
dists_order <- order(dists)[1]
min_dist <- dists[dists_order] # distance to closest
# the closest
sub_gen_tmp <- sub_gen[dists_order,]
sub_gen_tmp$min_dist <- min_dist
# weighted distance
sub_gen_tmp$He_harm_w <- weighted.mean(sub_gen$He_harm, 1/(dists^2)) # inverse of the distance
sub_gen_tmp$He_harm2_w <- weighted.mean(sub_gen$He_harm2, 1/(dists^2)) # inverse of the distance
sub_gen_tmp$He_harm3_w <- weighted.mean(sub_gen$He_harm3, 1/(dists^2)) # inverse of the distance
gen_data_v1_dist2_tmp <- rbind(gen_data_v1_dist2_tmp,
append(list(sub_gen_tmp),
list(sub_bio_area_j_tmp)) %>%
reduce(full_join, by="spp"))
}
gen_data_v1_dist2 <- rbind(gen_data_v1_dist2,
gen_data_v1_dist2_tmp)
}
gen_data_v1_dist2 <- gen_data_v1_dist2 %>% dplyr::filter(!is.na(He), !is.na(SHIFT))
write.csv(gen_data_v1_dist2,"Data/gen_data_final_m5.csv",row.names = FALSE)
# n range shifts
nrow(gen_data_v1_dist2)
## [1] 1881
# n genetic diversity measurements (after averaging by location)
dim(gen_div)
## [1] 7847 17
plot(gen_data_v1_dist2$He_harm,gen_data_v1_dist2$He_harm_w)
plot(gen_data_v1_dist2$He_harm2,gen_data_v1_dist2$He_harm2_w)
plot(gen_data_v1_dist2$He_harm3,gen_data_v1_dist2$He_harm3_w)
# N species
length(unique(gen_data_v1_dist2$spp))
## [1] 304
# Use only plus plus or minus minus
table(gen_data_v1_dist2$shift_vel_sign)
##
## negneg negpos posneg pospos
## 36 793 63 989
# N species
length(unique(gen_data_v1_dist2$spp))
## [1] 304
toplot_ <- gen_data_v1_dist2 %>%
group_by(Group) %>%
dplyr::summarise(Species = length(unique(spp)),
Obs = length(spp))
toplot_ <- reshape2::melt(toplot_)
## Using Group as id variables
ggplot(toplot_, aes(x = Group, y = value))+
geom_bar(stat="identity")+
geom_text(aes(y=value, label=value),
hjust = -.1, vjust=0.2, size=3,
position = position_dodge(0.9))+
theme_classic()+
coord_flip()+
facet_wrap(.~variable, scales = "free", ncol = 1)
toplot <- rbind(
data.frame(dataset = "Full",
lags = biov1$lag2,
shifts = biov1$SHIFT),
data.frame(dataset = "Subset",
lags = gen_data_v1_dist2$lag2,
shifts = gen_data_v1_dist2$SHIFT))
# lags
ggplot(toplot, aes(x = dataset, y=lags, fill = dataset, color = dataset))+
geom_boxplot(alpha = .5, outlier.shape = NA)+
scale_y_continuous(limits = quantile(toplot$lags, c(0.1, 0.9), na.rm = T))
## Warning: Removed 6280 rows containing non-finite values (`stat_boxplot()`).
tapply(toplot$lags, toplot$dataset, summary)
## $Full
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -145.1653 -0.4203 0.5907 0.3758 2.1868 13.1259
##
## $Subset
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -38.3203 0.0271 1.0197 0.9206 2.3096 13.1259
mod1 <- aov(lags~dataset, data = toplot)
summary(mod1)
## Df Sum Sq Mean Sq F value Pr(>F)
## dataset 1 525 525.5 24.45 7.66e-07 ***
## Residuals 31936 686320 21.5
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# shift
ggplot(toplot, aes(x = dataset, y=shifts, fill = dataset, color = dataset))+
geom_boxplot(alpha = .5, outlier.shape = NA)+
scale_y_continuous(limits = quantile(toplot$shifts, c(0.1, 0.9), na.rm = T))
## Warning: Removed 6376 rows containing non-finite values (`stat_boxplot()`).
tapply(toplot$shifts, toplot$dataset, summary)
## $Full
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -119.1696 -0.4356 0.2908 1.1661 2.4130 146.3000
##
## $Subset
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -119.1696 -0.3667 0.1429 0.9984 1.8912 39.3400
mod1 <- aov(shifts~dataset, data = toplot)
summary(mod1)
## Df Sum Sq Mean Sq F value Pr(>F)
## dataset 1 50 49.78 1.694 0.193
## Residuals 31936 938576 29.39
There are difference in mean values of lags, but no difference in mean shift values, between the full and subset datasets.
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = lagC))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3, aes(color = SHIFT_abs))+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## `geom_smooth()` using formula = 'y ~ x'
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
## -Inf
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## `geom_smooth()` using formula = 'y ~ x'
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
## -Inf
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## `geom_smooth()` using formula = 'y ~ x'
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
## -Inf
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = lagC, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = lagC, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## `geom_smooth()` using formula = 'y ~ x'
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
## -Inf
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## `geom_smooth()` using formula = 'y ~ x'
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
## -Inf
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## `geom_smooth()` using formula = 'y ~ x'
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
## -Inf
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = lagC2, color = SHIFT_abs))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = lagC2, color = SHIFT_abs))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("lagC (Velocity - Shift)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = log1p(sqrt(abs(SHIFT_cor))), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = log1p(sqrt(abs(SHIFT_cor))), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Param, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## `geom_smooth()` using formula = 'y ~ x'
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
## -Inf
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## `geom_smooth()` using formula = 'y ~ x'
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
## -Inf
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## `geom_smooth()` using formula = 'y ~ x'
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
## -Inf
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(.~Group, scales = "free", ncol = 3)
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm, y = (abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm, y = (abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm2, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm2, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "LAT") %>%
ggplot(aes(x = He_harm3, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Latitude")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
gen_data_v1_dist2 %>%
dplyr::filter(Type == "ELE") %>%
ggplot(aes(x = He_harm3, y = sqrt(abs(SHIFT_cor)), color = lag))+
ggtitle("Elevation")+
xlab("Genetic diversity (He harm logit)")+
ylab("Range shift (km/year)") +
scale_color_viridis_c(direction = -1)+
geom_point(alpha = .5, size = 3)+
theme_bw()+
geom_smooth(method = "lm", color = "black", size=1)+
stat_smooth(geom="line", method = "lm", alpha=0.3, se = FALSE,
aes(by=as.factor(Article_ID)), color = "black")+
facet_wrap(Param~Group, scales = "free")
## Warning in stat_smooth(geom = "line", method = "lm", alpha = 0.3, se = FALSE, :
## Ignoring unknown aesthetics: by
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'